'''simple dc model of the bsudaq analog stage
'''

import sys
import traceback
from math import log10
import numpy as np
import numpy.linalg as linalg
try:
    import matplotlib.pyplot as mpl
except:
    mpl = None

class AnalogStage(object):

    def __init__(self):

        self.Rd = 10e4
        self.Rs = 1e3
        self.Rf = 300e3

        self.Vs = 1e-3
        self.Vc = 5.0
        self.Aol = 1e6


def computeAOut(amp):
    '''variables are Vo, Vp, Vm.  equations are Gain, KVp, KVm.'''

    A = np.zeros((3,3))
    b = np.zeros((3,))
    A[0,0] = 1.0
    A[0,1] = -amp.Aol
    A[0,2] = amp.Aol
    b[0] = 0.0
    
    A[1,0] = 0.0
    A[1,1] = -2.0/amp.Rd - 1.0/amp.Rs
    A[1,2] = 1.0/amp.Rs
    b[1] = -amp.Vc/amp.Rd + amp.Vs/amp.Rs

    A[2,0] = 1.0/amp.Rf
    A[2,1] = 1.0/amp.Rs
    A[2,2] = -1.0/amp.Rf - 1.0/amp.Rs
    b[2] = amp.Vs/amp.Rs

    x = linalg.solve(A, b)

    return x


def dcGain(amp):

    Vs = 1.0e-3
    amp.Vs = Vs
    Vo1 = computeAOut(amp)[0]
    amp.Vs = - Vs
    Vo2 = computeAOut(amp)[0]
    return (Vo1 - Vo2) / (2.0 * Vs)



class Usage(Exception):
    def __init__(self, msg):
        self.msg = msg

class Error(Exception):
    def __init__(self, msg):
        self.msg = msg

def main(argv=None):
    import getopt

    if argv is None:
        argv = sys.argv

    options = "vh"
    verbose = 0

    try:
        try:
            opts, args = getopt.getopt(argv[1:], options, ["help"])
            for opt, val in opts:
                if opt == '-h':
                    print __doc__
                    return 0
                elif opt == "-v":
                    verbose += 1

        except getopt.error, msg:
            print >> sys.stderr, __doc__
            raise Usage(msg)

        x = []
        y = []
        amp = AnalogStage()
        print "        Rd   gain"
        for RdRange in [1.0, 10.0, 100.0, 1000.0, 10000.0, 100000.0]:
            for RdFac in [1.0, 2.0, 3.0, 5.0, 7.0, 9.0]:
                amp.Rd = RdFac * RdRange
                g = dcGain(amp)
                print >> sys.stderr, "%10.1e  %5.1f" % (amp.Rd, g)
                x.append(log10(amp.Rd))
                y.append((g/300.0))
        if mpl:
            mpl.plot(x, y)
            mpl.xlabel("log10(RD)")
            mpl.ylabel("G/G0")
            mpl.show()

    except Usage, err:
        print >>sys.stderr, err.msg
        print >>sys.stderr, "for help use --help"
        return 2

    except RuntimeError, err:
        print >>sys.stderr, err
        return 3

    except:
        traceback.print_exc()
        return 1

    finally:
        try:
            evdb.close()
        except:
            pass


if __name__ == "__main__":
    sys.exit(main())
    
